Mixing UTM and geographic axes annotations

This is a question that comes up regularly such that GMT has a gallery example (famous example 28) to show how to do it. But even with that example it's not a trivial thing to do.

However it's trivial from the Julia wrapper under the condition that the grid or image to be displayed is referenced internally. If it is not, grdedit can be used to assign a referencing system via its option proj. The other condition is that grid or image to be displayed is already in memory (so that the internal magicks can work). Given that, example 28 boils down to this (with a little bit less fancy details).

using GMT
G = gmtread("@Kilauea.utm.nc");
C = makecpt(cmap="copper", range=(0,1500));
imshow(G, cmap=C, shade=true, frame=(axes="WS", annot=true),
	coast=(shore=true, ocean=:lightblue, frame=(axes="EN", annot=true, grid=true)))

The above example, however, relied in having a UTM-referenced grid file already available. If you want to do the same with a global grid such as @earth_relief_03s, you can request a region in UTM coordinates and ask to project the fetched grid. This can be done on-the-fly by using the region, proj and the convert=true options in gmtread, like so:

using GMT
G = gmtread("@earth_relief_03s", region=[420000, 470000, 4510000, 4535000],
            proj="+proj=utm +zone=33", convert=true);
imshow(G, shade=true, frame=(axes="WS", annot=true),
	coast=(shore=true, ocean=:lightblue, frame=(axes="EN", annot=true, grid=true)))